Sticky obstacles to intramolecular energy flow 
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Vibrational energy flows unevenly in molecules, repeatedly going back and forth between trapping 
and roaming. We identify bottlenecks between diffusive and chaotic behavior, and describe generic 
mechanisms of these transitions, taking the carbonyl sulphide molecule OCS as a case study. The 
bottlenecks are found to be lower-dimensional tori; their bifurcations and unstable manifolds govern 
the transition mechanisms. 
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Chemical reactions usually proceed through a complex 
choreography of energy flow processes that deliver the 
needed vibrational energy to the reactive mode. The 
manner and time in which energy travels determine the 
outcome of the reaction and the properties of the prod- t,£ p 
ucts. The conventional wisdom concerning this funda- 
mental process is that vibrational energy travels very fast 
and well before a reaction takes place, distributes itself 
statistically among the modes of the molecule, assumed 
to resemble an ensemble of coupled oscillators. Reaction 
rate theories based on these assumptions - known collec- 
tively as statistical theories [l[ - have been vindicated in a 
number of chemical reactions. However, there is increas- 
ing evidence that the approach to equilibrium usually 
proceeds more slowly than predicted by statistical theo- 
ries Q - and it is also nonuniform, showing intriguing fits 
and starts. This anomalous diffusion is caused by vari- 
ety ofphase space structures, such as resonant islands or 
tori [3( that strongly slow down the trajectories passing 
nearby 0, 0| and therefore are said to be "sticky" [B| . To 
date, the theories so successfully applied in pioneering 
works @, 0, H, 0, E3] to lower-dimensional systems have 
not been extended beyond two degrees of freedom due to 
severe technical difficulties 11. 13. Ha. 
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FIG. 1: The generic behavior of chaotic trajectories in Hamil- 
tonian systems involves substantial fraction of intermittent 
behavior. The time-frequency analysis of a typical OCS tra- 
jectory (top panel) allows one to register the transition region 
(shaded band) and the frequencies £ of the regular motion, 
while the time series (lower panel) display the striking features 
of this abrupt change, t is time (in units of To = 0.063 x 10~ 12 
s) and £pj are the frequency ridges (in units of Tg -1 ) in the 
time- frequency decomposition [l4[ of Pi{t). 



The OCS molecule, an important player in greenhouse 
effect [151 ] . displays the slow and uneven relaxation to 
statistical equilibrium mentioned above, despite its three 
strongly coupled degrees of freedom. Models of OCS have 
served as a testbed for studying intramolecular dynam- 
ics in the chaotic regime [l6| and these classical findings 
have been confirmed in parallel quantal wave packet cal- 
culations [12| . In this Letter, we investigate vibrational 
energy flow in the OCS molecule using a Hamiltonian of 
the form [HI 

H = T(R U R 2 , a, P 1 ,P 2 ,P a ) + V(R 1 ,R 2 ,a), (1) 

where T is the standard kinetic energy of a rotationless 
triatomic molecule represented by two interatomic dis- 
tances i?i and R 2 , and a bending angle a (with their 
canonically conjugate momenta P\, P 2 and P a ). The 



potential V is fitted as 



V(R U R 2 , 
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K(ii i ) + Vr( J R 1)J R 2 ,i?3), (2) 



where R3 is the distance between O and S. The potential 
consists of Morse potentials Vi for each diatomic pair 
and an interaction potential V/ of the Sorbie-Murrell 
form 16]. A rich mixture of chaotic and regular dy- 
namics is observed at energies close to dissociation [l9j. 
The computations below were performed at 90% of the 
dissociation energy of the weaker bond. Trajectories in 
the vicinity of a specific periodic orbit with elliptic nor- 
mal stability are studied, focusing on their escape to the 
chaotic region, and identifying a generic mechanism of 
crossover from diffusion [20] to hyperbolicity and chaos. 

Figure Q] displays the time series of such a trajectory, 
initially close to the periodic orbit O a with period Tq (see 
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FIG. 2: Projections of the trajectory near a periodic orbit O a (with period To), analyzed in Fig.[T] The trajectory is represented 
in (Ri, R2) plots, broken down into segments, corresponding to the trapping stage (left panel) and chaotic stage (center panel). 
The bottleneck of transition from diffusion to hyperbolicity can be identified as a two-dimensional invariant torus (right panel.) 
The trajectory is sampled at fixed time intervals To/2. The orbit O a is shown as a solid curve in the center. 



Ref. [131). Figure^] shows salient features of capture (left 
panel) followed by an abrupt transition to chaos (cen- 
ter panel). An alternative view of the transition mecha- 
nism appears in Fig. [3] in terms of the Poincare section 
E : P a — 0, P a > 0, a < %. A boundary, marking the 
crossover from diffusion to hyperbolicity, can be identi- 
fied in terms of an invariant two-dimensional torus (right 
panel of Fig. [2]). Normal bifurcations of two-dimensional 
tori turn out to be the key ingredients in the transi- 
tion mechanism to hyperbolicity, as will be shown below. 
Generically, there are two stages in the dynamics of the 
trajectory. During the trapping stage (duration itrap)j 
the trajectory is close to (quasi-)periodic, following the 
unstable manifold of normally hyperbolic tori with very 
small positive Lyapunov exponent (in our case, A ~ 10~ 2 , 
thus explaining the observed trapping time itrap ~ A ). 
During the escape stage (duration t csc ] the shaded band 
in Fig. [T]and "tentacles" in Fig. [3]), the trajectory follows 
the unstable manifolds of the periodic orbit which is in 
3:5 resonance with O a (thick dots in Fig. [3]), with a signif- 
icantly larger Lyapunov exponent, leading to a fast tran- 
sition to the chaotic region of phase space (center panel of 
Fig. [21). These two time scales usually satisfy t tra p ^ ^osc- 
Observations of repeated trapping-escape-chaotic pro- 
cesses in relatively short trajectory segments (~ 10 3 Xo) 
provide evidence that these effects are prevalent. Dy- 
namical systems theory identifies structures with mini- 
mal hyperbolicity as key players in describing long-term 
features of the chaotic component of an attractor, and 
integral surfaces with small positive Lyapunov exponent 
are candidates for the "backbone" of the dynamics. Nor- 
mally hyperbolic invariant manifolds [2~fl ] have recently 
been implicated in the symbolic dynamics and phase 
space partition of higher-dimensional chaotic Hamilto- 
nian systems 22|, |23[, systems with small-dimensional 
saddles such as the "Crossed Fields" [2J| and the Re- 



stricted Three Body Problems [25 1 



Using a combination of trajectory diagnostic tools like 
Lyapunov maps (l9l . [26| . time- frequency analysis 14 1. 
and methods from the theory of dynamical systems like 



periodic and quasiperiodic orbit computations 27|, [28(, 
we relate the phenomenon of trapping to invariant struc- 
tures in phase space and to lower-dimensional invariant 
tori (with a relation to their normal stability properties) 
in particular. It is commonly assumed that in "typical" 
Hamiltonian systems with a large number of degrees of 
freedom TV, the relative measure of A^-dimensional invari- 
ant tori (N local integrals) is either zero or one [29( . The 
implication is that chaotic systems with large N approach 
conditions of the stochastic ansatz, and hence, the trap- 
ping phenomenon described above is insignificant. On 
the other hand, it has been established recently that high 
order resonances form robust islands of secondary struc- 
tures with positive measure [30J. 

In order to identify bottlenecks of transition from diffu- 
sion to chaos, we monitor the progress of invariant phase 
space structures along the transition channel using ro- 
tation numbers. The results are summarized in Fig. 31 
which is central to understanding this transition. In a 
trapping region around the elliptic periodic orbit O a (left 
panel of Fig. [21), the rotation numbers are obtained from 
the frequency map analysis (2(3] on the surface of sec- 
tion. It can be characterized by a single w tra p ~ 0.60556, 
implying that a two-dimensional torus is the relevant in- 
variant structure in the trapping process. Having com- 
puted a family of two-dimensional tori, parametrized by 
rotation numbers uj, it is evident that w tra p places the 
torus on the hyperbolic branch of the bifurcation dia- 
gram represented in Fig. [U This implies that the es- 
cape is mediated by manifolds of a torus with hyperbolic 
normal stability. The duration of the trapping stage is 
approximately 150 returns on X, and is consistent with 
the maximal Lyapunov exponent A < 0.05. Processes 
associated with the escape from the trapping region can 
be better understood by analyzing the tangent space of 
the elliptic periodic orbit O a that locally has the struc- 
ture of a direct product (center + center) T x L x T x I 2 , 
with the periodic orbit at the origin. The elements of the 
two intervals I, C M are rotation numbers uji, which are 
not unique in general: The choice is fixed by requiring 
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where fi is a measure of the torus and 



w° are stability angles of the elliptic periodic orbit O a 
(lo° = 0.24500633 and lo° = 0.37046872). The Poincare 
map induces rotations on T, r ui x 1 x r U2 x 1, where r u 
is a rotation on T with the rotation number lo. Partial 
(or complete) resonances are determined by one (or two) 
resonance conditions ruox + uilo^ + k = 0, where (n, m, k) 
are integers such that \n \ + \ m\ + \ k\ > 0. The most strik- 
ing trapping effects are observed for partial resonances of 
the type Txli X {0} x {0}, and {0} x {0} x T x I 2 . Choos- 
ing either of the two situations, a resonance channel has 
been constructed by finding the two-dimensional invari- 
ant tori for Wj € Ij . In order to find these tori we consider 
the Poincare map J-^: : X i— > S. Tori may have hyper- 
bolic normal linear stability, therefore a search for them 
cannot rely on methods exploiting "stickiness" proper- 
ties. The sections of two-dimensional invariant tori are 
one-dimensional closed curves (called hereafter "loops"). 
We consider loops as discretizations of 7 : T 1— > S (with 
periodic boundary condition j(s) = j(s + 1)) and re- 
quire that the Poincare map T-z , restricted to the loop is 
equivalent to a rigid rotation r u . This translates into an 
invariance condition: 



-?s(7(s)) = j(s + u)). 



(3) 



Equation ([3]) is solved using damped Newton iterations 
for the Fourier coefficients of 7(5). The linear stability 
properties of the loop are determined by (A, -0), solutions 
of the generalized eigenvalue problem: 



Dfx(s)Ms) = AiP{s + lu). 



(4) 



Equation Q has a one-dimensional kernel, which we 
eliminate using singular value decomposition. The initial 
data for the Newton iterations 70 (s) and to were obtained 
using one of the following two methods : The first method 
uses the trapping region of the trajectory (see Fig.]]}. We 
estimate lo using Fourier- like methods [201 ] , and truncate 
the continued fraction expansion of lo = [ai, 0,2, ■ ■ •] be- 
fore the first large Oj so that loq — P/Q. Then we take 
sequences of trapping region data every Q iterations and 
combine them to obtain 70 (s). A refined value of to can 
be estimated by minimizing | J-^ o 70 — 70 o r u \ . The second 
method combines continuation in lo with the direct prod- 
uct structure in the neighborhood of the periodic orbit. 
The surface of section derivative DTy, at the periodic 
orbit has two pairs of complex eigenvalues exp[±ta>°], 
i = 1,2. The eigenvectors define mutually skew orthogo- 
nal symplectic vector spaces Vt = R 2 . It is assumed here 
that the linear approximation is effective in the neigh- 
borhood of the periodic orbit. 

The set of two-dimensional tori is found to be dis- 
continuous at the gaps in Fig. |4] due to complete reso- 
nances (periodic orbits) and secondary invariant struc- 
tures. Normal stability is typically elliptic for small 
\lo — lo®\. We identify the two-dimensional invariant torus 
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value 


Cont. frac. 


LO° A 


0.240711317575 


[4, 6, 2, 11, 5, 5...] 


LO% 


0.215852976389 


[4, 1, 1, 1, 2, 1, 1, 1, 1, 2, 1, 1...] 


c 

LO C 


0.608654398762 


[1, 1, 1, 1,4, 45, 1, 1,1,1...] 


c 

LO D 


0.605804087926 


[1, 1, 1, 1, 6, 3, 2, 2, 1...] 



TABLE I: Rotation numbers of the two-dimensional invariant 
tori at the bifurcation points A, B, C, D shown in Fig. 4. 
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FIG. 3: Poincare section of the trajectory near a periodic 
orbit O a , analyzed in Figs. [T] and [5] The bottleneck (a two- 
dimensional torus) is a loop (blue) at the bifurcation point 
("D" in Fig. [1J. The trajectory is trapped in the vicinity 
of a loop (which is clearly seen from the inset). The escape 
stage is shown as two "tentacles," which extend along the 
unstable manifolds of a resonant periodic orbit (the five red 
dots around the center). 



at the period doubling bifurcation point as a bottleneck 
of a given resonance channel. The rationale follows from 
the theory of dynamical systems: Beyond the bifurca- 
tion point at lo = lo c , the normal stability changes to 
hyperbolic. This change affects trajectories passing by 
its neighborhood. One recurrent observation is that the 
continued fraction expansion of bifurcation rotation num- 
bers has a tail composed of small integers (see Tab. HJ. 
This feature is reminiscent of the observation that the 
continued fraction expansion of the frequency of the last 
invariant torus in generic Hamiltonian systems with two 
degrees of freedom is noble (with a tail of ones) in many 
situations [6j. 

The reliability of the numerical solution can be tested 
by examining its Floquet multipliers, given by Eq. ([?]). 
An exact solution consists of a set of complex numbers 
with up to three different absolute values: 1, A, 1/A. Sig- 
nificant variation from these values signals an unreliable 
solution. 

In conclusion, our findings indicate that trapping and 
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FIG. 4: Fine structure of invariant tori, scanned along the 
transition channel. The plot shows how Lyapunov exponents 
depend on the rotation number u. The points of frequency 
halving bifurcations ("A"- "D") can be interpreted as bottle- 
necks of transition from diffusion to hyperbolicity. Red dots: 
family of loops arising from the periodic orbit O a - Black dots: 
frequency halved loop, emerging at the bifurcation point "A". 
Insets display (7?i,Pi) projections of loops near the bifurca- 
tion point "A". Red: loop with elliptic normal stability and 
uj — uji ~ 0.24067. Black: loop with hyperbolic normal sta- 
bility and uj = uj 2 = (uji + l)/2 « 0.62033. 



escape are mediated by the same sequence of events, 
and an approximate boundary, which separates trapped 
and chaotic behavior, can be found in analogy with 
the boundaries that separate reactants from products in 
Transition State Theory [l| , where sharply defined phase 

M 



space structures 




play this role. 



In a broader context, our work forms yet another stim- 
ulus to reconsider the relevance of local integrals and par- 
tial resonances in realistic, chaotic Hamiltonian systems 
with many degrees of freedom. Here, we have explained 
a paradoxical situation, namely that integral surfaces 
with positive Lyapunov exponents (i.e., not "sticky") 
can trap chaotic trajectories. Widespread observations of 
repeated trapping-escape-chaotic processes in short tra- 
jectory segments provide evidence that these effects are 
generic and occurring frequently in many settings rang- 
ing from plasmas to celestial mechanics. 
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